Weak Disorder in Fibonacci Sequences 
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We study how weak disorder affects the growth of the Fibonacci series. We introduce a family 
of stochastic sequences that grow by the normal Fibonacci recursion with probability 1 — e, but 
follow a different recursion rule with a small probability e. We focus on the weak disorder limit 
and obtain the Lyapunov exponent, that characterizes the typical growth of the sequence elements, 
using perturbation theory. The limiting distribution for the ratio of consecutive sequence elements is 
obtained as well. A number of variations to the basic Fibonacci recursion including shift, doubling, 
and copying are considered. 
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The Fibonacci integer sequence {1,1,2,3,5,8,...} has 
been studied extensively in number theory, applied math- 
ematics, physics, computer science, and biology Q, Q. 
Fibonacci numbers are ubiquitous in nature: they gov- 
ern branching in trees, spiral patterns in shells, and the 
arrangement of seeds in sunflowers 0, 0, 01 ■ 

The Fibonacci sequence, defined via the recursion re- 
lation 



F n +x — F n + F n -i 



(1) 



with Fq = and Fx = 1, is deterministic. However, many 
patterns in nature are not perfect. For example, spiral 
patterns in sunflowers, where Fibonacci numbers as high 
as 144 are observed, may very well be disordered. An em- 
pirical study of sunflowers observes the normal sequence 
{1, 1, 2, 3,5,.. .} with a frequency of 95%, but altered se- 
quences such as {2, 3, 5,7,.. .} and {1, 3,4,.. .} are also 
observed with a small frequency . 

Motivated by this empirical observation, we study dis- 
order in Fibonacci sequences. Specifically, we introduce 
the following stochastic sequence. We assume that the 
normal Fibonacci rule is followed most of the time, 
but that with a small probability, e< 1, an "erroneous" 
recursion x n +x = x n + x n -i, involving an index shift, is 
followed. The stochastic recursion rule is therefore 



<-n+l 




x n -x prob 1 - e; 
X n -% prob e. 



(2) 



The initial elements are always xq = 0, and Xx = x% = 1. 

Let us first recall a few useful facts on the Fibonacci 
sequence, corresponding to the limiting case e = 0. The 
series elements are given by 



F n = 



A" - (-A)- 
A + A" 1 



(3) 
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with the golden ratio A = ■ ^he ser i es elements 

grow exponentially with n, F n ~ A". Substituting this 
form into the recursion £[J, the golden ratio satisfies 
A 2 = A + 1 , and this allows to express polynomials of 
arbitrary degree in A as linear functions of A. More- 
over, the ratio between two successive series elements, 
f n = x n /x n -i, approaches the golden ratio r n — > A, as 
n — > oo. 

Our goal is to elucidate the typical growth of the se- 
quence elements 



(4) 



with (3 = /3(e) the Lyapunov exponent. For exam- 
ple, for the random Fibonacci S6I"1GS X ji — — 1 — L Xyi — 2 
where addition and subtraction are chosen with equal 
robabilities, the Lyapunov exponent is j3 rs 0.123975 
H, 0, 0, [n|. In our case, however, the recursion 
rules 10) represent a gentle departure from the original 
Fibonacci rule Q and thus, we expect a small change in 
the Lyapunov exponent. We focus on the weak disorder 
limit, e — ► 0, and use perturbation theory to show that 
the Lyapunov exponent varies linearly with the disorder 
strength 



/3(e) =(30 + 13x6 + 



(5) 



with Pq = In A. 

In general, the average behavior (x n ) can be obtained 
analytically. From the sequence definition @, the aver- 
age satisfies the recursion relation 



(x n+1 ) = (x n ) + (1 - e)(a;„_i) + e(x n - 2 ) 



(6) 



with (xo) = 0, and (xx) = (^2) = 1- This linear rela- 
tion implies the exponential growth (x n ) ~ fj, n with the 
growth factor /x being the largest root of the third order 
polynomial 



/! 2 + (l-e)/i + e. 



(7) 
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Differentiating this equation with respect to e and setting 
e = 0, we find (i/Lt/(ie| e , and thus, for small e we have 



2 



jLt(e) = A — x+^-e- To compare with the growth of the 
typical sequence Q), it is useful to write (x n ) ~ e 7n with 
7 = ln/x. To first order in the disorder strength e, 



7( e ) = 7o + 7i £ 



(8) 



l-A 



with 70 = 0q and 71 = A(A+2 ) . 

To address the typical behavior, we introduce the ra- 
tio between two successive elements in the sequence, 
r n = x n /x n -\. The random recursion rule J3J implies 
that this ratio satisfies the random map, 



r n +i 




1 



(9) 



When there is no disorder, e = 0, the ratio approaches 
the golden number, r n — > A as n — ► 00. Thus as the num- 
ber of iterations of the normal map (0 grows indefinitely, 
the distribution of the ratio approaches a delta function 
centered at the golden ratio, P(r) — ► 5(r — A). 

Generally, when e > 0, the distribution P(r) has a 
richer structure, as shown below. The Lyapunov expo- 
nent can be conveniently expressed in terms of P(r). In- 
deed, each sequence element is given by the product 



n 



(10) 



With the exponential growth , the Lyapunov exponent 
simply equals the expected value of the logarithm of the 
ratio, 



0- 



[bxr] 



dr P(r) lnr. 



(11) 



At weak disorder, with a small probability e, an er- 
ror occurs. That is, the map r„+i = 1 + \/{r n r n -i) is 
implemented. As long as no errors occur, the ratio will 
essentially be equal to A. Then, when an error occurs, 
the ratio reduces to 1 + A~ 2 . Since the expected num- 
ber of iterations before another error occurs, e _1 , is very 
large, the ratio again quickly approaches A. This cycle 
continues ad-infinitum. 

To characterize this process, we should understand how 
an error evolves under the random map @. Thus, we 
consider the following scenario: (1) Initially, the ratio 
equals the golden number po = A, (2) An error occurs 
at the very first step, and (3) no further errors occur. 
Let p n be the value of the ratio after n iterations. Then 
pi = 1 + A~ 2 , and using the relation A 2 = A + 1 we have 
pi = jt^- At further iterations, the ratio follows the 
normal map p n +\ = 1 + 1/ Pn and therefore, 



Pi = 



2 + A 
1 + A 



2A 



3A 



P2 = 



P3 



3 + 2A 



By induction, at the n iteration, the ratio can be ex- 
pressed in terms of the Fibonacci numbers 



Pn+2 + Pn+lA 



+ 1 



(12) 



This series alternates around A: p2n+i < A while pi n > A, 
but both the odd and the even sub-series quickly converge 
to golden ratio, p n — > A as n — > 00. 

This analysis characterizes how a single error affects 
the ratio. To first order in the disorder strength e, 
the probability that the value p n is observed equals 
e(l — e)" _1 , reflecting the probability that one error is 
made and then, no errors are made in the following n—1 
iterations. Our first result is the distribution P(r) for 
the ratio to have the value r 



P(r) - 5 (r-pn), 



(13) 



71=1 



in the weak disorder limit e — > 0. 

The calculation of the first order correction to the Lya- 
punov exponent JSJ) is now straightforward. Substituting 
the leading behavior in the weak disorder limit l|13|) into 
the general formula (|llf> for the Lyapunov exponent we 
obtain 



A 



e)"- 1 lnp„ 



(14) 



The sum is evaluated as follows 



A^e^l-er^lnA + e^l-er- 1 !!!^. 



Performing the summation in the first term, we verify 
that Po = In A. Keeping only the terms proportional to 
e in the second sum gives the leading correction in the 
perturbation expansion of the Lyapunov exponent J5j), 



E ln - 



n=l 



A 



(15) 



To perform this summation, we substitute the expression 
(|12|l . and replace the upper limit with a large but finite 
cutoff N, and then evaluate the N — > 00 limit as follows 



01 



,. , 1 Fn+2 

hm m T-T7 



Fjy + i\ 



1 + A 



In- 



2A 2 



(A + 1)(A + A-!) 



where in the second line we used Eq. ©. Using the 
equality A 2 = A + 1 , we arrive at our second main result 



01 



In 



2A 
A"+2 



(16) 



This correction is very close, but not identical, to that 
corresponding to the average behavior (JSJ. As the dif- 
ference 0i — 71 w —0.00599897 is negative, the typical 
growth is slower than the average growth. This mani- 
fests the multiscaling behavior that has been reported in 
other stochastic sequences ^2 ■ Generally, there is a mul- 
tiscaling spectrum £ m that characterizes the growth of 



3 
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FIG. 1: The weak disorder limit. Shown is the Lyapunov ex- 
ponent (3, versus the disorder strength e. Also shown for ref- 
erence are the parameter 7 characterizing the average growth 
and the perturbation theory result Jj)J with 1161 . 



the mth moments, (a;™) 1 /™ ~ exp(n£ m ). However, there 
is no obvious relation between the Lyapunov exponent /3 
and the multiscaling spectrum £ m , e.g., (3 ^ 7 = £i- 

We performed Monte Carlo simulations to verify these 
theoretical predictions. In the simulations, we followed 
the stochastic evolution of the variable r. This approach 
is advantageous for computation because the ratios are 
bounded, in contrast with the explosive growth in the se- 
quence elements. The results presented here correspond 
to a single Monte Carlo run with 10 9 iterations. 

There is a distinct but subtle difference between the 
typical and the average growth as characterized by j3 and 
7, respectively (Fig. 1). The two coincide in the limit- 
ing cases e = and e = 1 p|, and the discrepancy is 
maximal, a mere 0.2%, at the midpoint e = 1/2. 

The numerical simulations show unambiguously that 
as the number of iterations grows indefinitely, the ra- 
tio distribution approaches a stationary distribution 
P{r). On Fig. 2, we display the cumulative distribution 
G(r) = f^dr'P(r'). 

The stationary distribution has a compact support, 
''mm < r < r max . Indeed, the definition of the map @ 
implies the obvious bounds r m - m > 1 and r max < 2. The 
values r m ; n = (l + y3)/2 and r max = \/3, consistent with 
the numerical simulations results, are obtained from the 
following relations 



1 



1 



l + r n 



(17a) 
(17b) 



The first relation (|17afl follows from the normal Fibonacci 
recurrence r„+i = 1 + l/r„. The second relation l|17b(l 
follows from the altered recurrence r„+i = l + l/(r„r„_i) 
combined with max(r„r„_i) = l+max(r„_i) = l + r max , 
that follows from the normal recursion r n+ \ = l + l/r„_i. 
The distribution P{r) consists of a set of delta func- 




.3 1.4 1.5 1.6 1.7 1.8 
r 

FIG. 2: The cumulative distribution G(r) versus the ratio r. 



tions, and therefore, the cumulative distribution G(r) 
has a devil's staircase structure with infinitely many 
gaps. Generally, there is a large gap in the interval 
l/2<r<l + l/-\/3- This gap arises since the map 



@ transforms (r 



mini ''max 



) into the union of two subin- 



tervals, (r min ,3/2) and (1 + l/r raax , r max ). The bounding 
point 3/2 is obtained using reasoning similar to that used 
in the previous paragraph. Restricting the map to the 
above subintervals one finds that they are transformed 
into the union of four smaller subintervals, etc. Hence the 
support of the invariant distribution P(r) is a Cantor-like 
fractal set and the cumulative distribution therefore has 
a devil's staircase structure with an uncountable number 
of singularities (Fig. EJ) . 

Thus far, we have addressed a specific modification 
of the Fibonacci recurrence, namely, the one involving 
the index shift x n +i = x n + x„_2- But there are of 
course several other, equally natural, modifications of the 
basic recursion rule. For example, one may simply copy 
the last element x n+ \ — x n or alternatively, double it 
x n+ \ = x n + x n . These two models are analyzed along 
the same lines. For simplicity, we address the latter case 
where the recursion relation is 




prob 1 — e, 
prob e; 



(18) 



with the initial elements x\ 
random map is 



Xl 



1. The corresponding 



'n+l = 




prob 
prob 



1 - 
e. 



e; 



(19) 



In contrast with 10, when an error occurs, the ratio 
r = 2 is independent of the previous element. Thus, er- 
ror events effectively reset the process anew. As a result, 
this stochastic process is analytically tractable. 

To characterize how an error propagates, we start with 
pi = 2 and use the recursion p n +i = 1 + 1/p« to obtain 
the first few terms P2 = 3/2, p% = 5/3, and p^ — 8/5. In 



4 



general, 



F, 



n+2 



F 



n+1 



The ratio attains this value, r — p n , when an error is 
followed by n — 1 normal recursion steps, and this oc- 
curs with probability e(l — e)™ -1 . Thus, the probability 
distribution of the ratio is 



P(r) 



^(1 



(21) 



with p n given by (|20[1 . In contrast with the limiting dis- 
tribution l|13l) . this result is now exact, because the his- 
tory prior to the most recent error event is irrelevant. 
The distribution now has a countable set of singularities 
located at the ratios p n of successive Fibonacci numbers. 
These singularities "bunch" near the golden ratio A. 

Substituting the probability distribution <|21[) into the 
Lyapunov formula yields 



/9 = e£(l 



\n-l 



hi 



F 



n+2 



Again, the typical growth is slower than the aver- 
age growth, as for example, /?(l/2) « 0.571357 while 
7(1/2) w 0.577049 14]. The exact expression |g2J can 
be, in principle, expanded as a power series in the disor- 
der strength e, viz. = J2 n >o /3n£ n with (3 n characteriz- 
ing the effect of n errors. OFcourse, (3q = hi A. The low- 
est order correction, which can be obtained either from 
Eq. or from Eq. |JTSJ, is given by 



In- 



2A + 1 
A + 2 



One can also extract the next correction from Eq. 122|l : 
the result is fa = E m >o ln t 1 + (-l) m A- 2ro - 6 ] . 

(20) In summary, we introduced a class of random Fi- 

bonacci sequences where with a fixed probability the clas- 
sic rule is followed, but otherwise, an alternate recursion 
occurs. We analyzed the weak disorder limit and ob- 
tained the limiting distribution for the ratio of consec- 
utive sequence elements as well as the Lyapunov expo- 
nent. We found that the typical growth is slower than 
the average growth. We also showed that the cumulative 
distribution of the ratio of consecutive elements has a 
devil's staircase structure. An exact solution for particu- 
larly simple alterations of the recursion rule was obtained 
as well. 

The above results raise a number of questions: Can 
the ratio distribution and the Lyapunov exponent be ob- 
tained analytically in general? What are the locations of 
the singularities underlying the distribution of the ratio? 
What is the probability that a given integer belongs to 
the random Fibonacci sequence? 

We have focused on the average and the typical se- 

(22) quence growth, but further information is encoded in 
fluctuations with respect to the typical behavior. Related 
studies on disordered systems suggest that such fluctu- 
ations should obey Gaussian statistics and our pre- 
liminary numerical simulations support this. The corre- 
sponding variance may be calculated using perturbation 
theory in the weak disorder limit. 
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